Paleopteran molecular clock: Time drift and recent acceleration

Abstract Applying BEAST v1.10.4, we constructed a Bayesian Inference tree comprising 322 taxa, primarily representing Paleoptera (Odonata and Ephemeroptera; Pterygota), Zygentoma and Archaeognatha (Apterygota; paraphyly), and Neoptera (Plecoptera; Pterygota), based on a 2685 bp sequence dataset. Our analyses revealed that robust dating required the incorporation of both Quaternary and pre‐Quaternary dates. To achieve this, our dating incorporated a 1.55 Ma (Quaternary) geological event (the formation of the Ryukyu Islands) and a set of chronologically well‐founded fossil dates, spanning from up to 400 Ma (Devonian) for the stem Archaeognatha, 320 Ma (Carboniferous) for the crown of Paleoptera, 300 Ma (Carboniferous) for the crown Ephemeroptera, and 280 Ma (Permian) for the crown Odonata, down to 1.76 Ma (Quaternary) for Calopteryx japonica, encompassing a total of 22 calibration points (events: 6, fossils: 16; Quaternary: 7, pre‐Quaternary: 15). The resulting dated tree aligns with previous research, albeit with some dates being overestimated. This overestimation was mainly due to the lack of Quaternary calibration and the exclusive dependence on pre‐Quaternary calibration, though the application of maximum age constraints also played a role. Our minimum age dating demonstrates that the molecular clock did not uniformly progress, rendering rate dating an inapplicable approach. We observed that the base substitution rate is time‐dependent, with an exponential increase evident from around 20 Ma (Miocene) to the present time, exceeding an order of magnitude. The extensive radiation and speciation of Insecta and Paleoptera, potentially resulting from the severe climatic changes associated with the Quaternary, including the commencement of glacial and interglacial cycles, may have significantly contributed to this increase in base substitution rates. Additionally, we identified a potential peak in base substitution rates during the Carboniferous period, around 320 million years ago, possibly corresponding to the Late Paleozoic Ice Age.

(Quaternary) for Calopteryx japonica, encompassing a total of 22 calibration points (events: 6, fossils: 16; Quaternary: 7, pre-Quaternary: 15).The resulting dated tree aligns with previous research, albeit with some dates being overestimated.This overestimation was mainly due to the lack of Quaternary calibration and the exclusive dependence on pre-Quaternary calibration, though the application of maximum age constraints also played a role.Our minimum age dating demonstrates that the molecular clock did not uniformly progress, rendering rate dating an inapplicable approach.
We observed that the base substitution rate is time-dependent, with an exponential increase evident from around 20 Ma (Miocene) to the present time, exceeding an order of magnitude.The extensive radiation and speciation of Insecta and Paleoptera, potentially resulting from the severe climatic changes associated with the Quaternary, including the commencement of glacial and interglacial cycles, may have significantly contributed to this increase in base substitution rates.Additionally, we identified a potential peak in base substitution rates during the Carboniferous period, around 320 million years ago, possibly corresponding to the Late Paleozoic Ice Age.

K E Y W O R D S
atmospheric CO 2 , BEAST v1.10.4,carboniferous, exponentially increased base substitution rate, fossil and geological event calibration, glacier period, Odonata, quaternary, tMRCA; quaternary and pre quaternary calibration

| INTRODUC TI ON
The recent trend in phylogenetic studies may be leading to larger tree and genome sizes, potentially resulting in more precise and detailed phylogenetic trees (Table 1).However, while increased sizes can enhance the robustness of topology, they do not guarantee strict dating accuracy.In other words, tree and genome sizes alone cannot calibrate the phylogenetic tree and are therefore inadequate for dating purposes (cf., Cicconardi et al., 2023).We perform node dating robustly and present a well-established Bayesian Inference (BI) dated tree of Paleoptera (Figure 1), examining its implications, with a specific focus on the methodological approach of correlating base substitution rates with time, similar to Osozawa (2023).BEAST v1.10.4 includes a feature that allows for the visualization of base substitution rates and ages on each node (c.f., Cicconardi et al., 2023), facilitating the exploration of time dependencies (Figure 1 inset).
In the earlier phase, a molecular clock model could estimate the age of a tree by assuming a relatively constant base substitution rate.
For instance, a rate of 0.0115 substitutions per site per million years (s/s/my) was established for adaptive radiated Heliconius butterflies and is considered standard (Brower, 1994).Papadopoulou et al.
Especially concerning sister-related species, it is crucial to note that in the maximum likelihood (ML) tree, each branch length, representing the number of substitutions per site, is distinct (e.g., see figure 2 in Osozawa, Sato, et al., 2017).This variability indicates that the base substitution rate progresses differently, rendering the application of the molecular clock hypothesis, or rate dating less stringent.One potential solution for molecular dating is to consider a relaxed clock model instead of a strict clock model (Zhang & Drummond, 2020), and dating applications offer the option to select the relaxed clock.
Molecular evolutionary rates, as seen in primate mitochondrial genes, have been observed to display both high short-term (1-2 million years) mutation rates and low long-term substitution rates, with the latter possibly attributed to purifying selection (Ho et al., 2005(Ho et al., , 2011;;BEAST v.1.3;Drummond & Rambaut, 2003).This timedependent nature of the molecular clock can influence a span of up to 1 million years (Papadopoulo et al., 2010).Although these studies have identified an increased intra-specific genealogical rate during the Quaternary period, they received limited attention, with the exception of the works of Osozawa and Wakabayashi (2022;latest version 2023;BEAST v1.10.4;Suchard et al., 2018) and Osozawa (2023;BEAST v1.10.4).
In our pursuit of revisiting the time-dependent phenomenon of increased base substitution rates during the Quaternary period (Ho et al., 2005;Papadopoulo et al., 2010), we have leveraged the Quaternary calibration date of 1.55 ± 0.15 Ma provided by Coeliccia damselflies endemic to the Ryukyu Islands (Osozawa et al., 2012;Osozawa, Sato, et al., 2017).We applied the same calibration method to endemic Platypleura in a global cicada BI tree (Osozawa, Shiyake, et al., 2017;Osozawa & Wakabayashi, 2022, updated in 2023).In TA B L E 1 Dated tree compilation.this paper, we extend this approach by incorporating an additional five Ryukyu endemic Odonata species, all of which make use of the same calibration date of 1.55 ± 0.15 Ma (Figure 1; Q1-Q6).
Osozawa and Wakabayashi (2022) focused on Cicadoidea, including Tettigarctidae (Hemiptera; Neoptera).In this paper, our focus is on Paleoptera, specifically Ephemeroptera and the primary subject, Odonata.In the Insecta classification, we find the paraphyletic group referred to as "Apterygota" and the clade Pterygota, and as outgroups, we have incorporated the "apterygotan" Archaeognatha (jumping bristletail) and Zygentoma (silverfish).
Pterygota includes Palaeoptera, which lacks the ability to fold the wings back over the abdomen, in contrast to Neoptera (Ishiwata et al., 2011).Additionally, we have included Plecoptera (stonefly) as one of the primitive groups of Neoptera, based on research by Misof et al. (2014) and Montagna et al. (2019).

| Materials
The total number of taxa in our analysis is 327, which is comparable to that of recent related studies listed in Table 1.

| DNA extraction and polymerase chain reaction amplification, sequence alignment
For the present study, we obtained a total of 58 sequence datasets,  Osozawa, Sato, et al. (2017), which also included a total of 30 datasets available in the present analyses, in addition to the previously published data in Osozawa and Wakabayashi (2015).
The COI-3P and COII-5P gene data overlap with the tRNA-Leu region, as described in Osozawa, Sato, et al. (2017), and note that data for the corresponding region are rarely available in GenBank/ DDBJ.However, COI-5P data, amplified using the LCO-HCO universal primers (Folmer et al., 1994;Osozawa, Takáhashi, et al., 2017), are frequently found in GenBank/DDBJ but are not compatible with our COI-3P sequence.
In BEAST v1.10.4, the analysis of combined genes, including mitochondrial COI, COII, and 16S rRNA, along with the nuclear 28S rRNA, can be performed simply by configuring partitions for these genes.Concatenation of these genes and re-partitioning are not necessary (Osozawa, 2023).
Our selected mitochondrial genes exhibit a high base substitution rate and higher resolution compared to nuclear genes, making them suitable for studying very old ages without experiencing mutation saturation, as discussed in Osozawa et al. (2012); Osozawa, Sato, et al. (2017), Osozawa and Wakabayashi (2022, updated in 2023), and Osozawa (2023).This is in contrast to the challenges posed by paralogs observed in nuclear genes (Gabaldón & Koonin, 2013).
Notably, Ho et al. (2005) found it improbable that the apparent decline in rates over time could be attributed to mutational saturation, and we propose the absence of mutation saturation (Osozawa & Wakabayashi, 2022, updated in 2023) and Osozawa (2023).To further explore this, we examined the relationship between pairwise distance and the number of transitions or transversions for each gene (Figure 2), applying the MEGA11 function (Tamura et al., 2021).The rapid subsidence required for each island's creation led to their simultaneous isolation from the Chinese mainland and from one another, primarily due to the formation of the Okinawa Trough and other significant seaways such as the Tsushima and Taiwan straits.

|
These straits, along with some minor ones in specific instances, We concur with their approach, which emphasizes reviewing the published geological age and the stratigraphic range of the fossil, and aligning the data with the most up-to-date standard global geochronology.Furthermore, they stress the importance of evaluating the quality of stratigraphic information for each fossil.Notably, they made exceptions for certain amber deposits, excluding most of them from their analysis.For instance, the Baltic amber was dated to the Lutetian period using Ar-Ar dating at 44.3 ± 0.4 Ma (Wappler, 2005) 2).However, it is essential to acknowledge that there is some degree of uncertainty surrounding fossil ages (not fossil identifications).Therefore, we conducted a reevaluation of the geological constraints linked to the fossil locations.We selected several pre-Quaternary calibration dates, taking into account their geological reliability, which included factors such as whether the fossil ages were constrained by modern and more precise radioisotopic dating or by stratigraphical correlation.
Table 3 provides a summary of our fossil calibration, with further details available in Appendix S2: fossil calibration.

| Picking the target tree with rate data
The target tree, fully calibrated by both fossil and geological event dates, and encompassing both pre-Quaternary and Quaternary dates, is presented in Figure 1.
We are revisiting a time-dependent phenomenon, specifically the exponential increase in phylogenetic base substitution rates F I G U R E 1 Bayesian Inference tree calibrated using both older fossil dates and the Quaternary geological event date with the aid of BEAST v1.10.4.Calibration points are represented by star marks and are detailed in Table 3.The calibration points labeled from A to K, and J (Quaternary) and X are based on fossil calibrations, while those marked Q1 to Q6 correspond to Quaternary geological event calibrations, specifically at 1.55 ± 0.15 Ma.For instance, Calibration Point I denotes the crown node of a common ancestor of the ingroup species of Libellulidae, but the stem node can be specified if necessary.In BEAUti, it is possible to specify a monophylum (clade) for Libellulidae ingroup species.It is important to note that the ingroup of geological calibration point Q6, for example, includes both Chinese Platycnemididae species and Ryukyu-Taiwan endemic species, as shown in Figure 3. Notably, some sequence data were incomplete or not whole mitochondrial.Marked as WM are sequences with complete whole mitochondrial data, including both COI-3P and COII-5P, as exemplified by data from Wang et al. (2021).In cases where the data are not entirely whole mitochondrial, available COI-3P, COII-5P, or both are indicated (marked in).For instance, Kim et al. (2014) included 43 taxa from Korea, and while COI-3P data are available, COII-5P data were not presented (marked #).In addition, mitochondrial 16Sr RNA data for Japanese taxa were complete, as submitted by Ozono et al. (2012).Sequence data retrieved from Carle et al. (2015) and Ware et al. (2014Ware et al. ( , 2017, concatenated gene) , concatenated gene) were incomplete and largely unavailable for the present BEAST v1.10.4 analyses (combined gene).Specifically, COI-3P and/or COII-5P data were lacking for species including Neopetalia punctata (Neopetaliidae, South America), Hemiphlebia mirabilis (Hemiphlebiidae; Australia), Heliocharis amazona (Dicteriadidae, South America), and Oligoaeschna (Sarasaeschna) pryeri (Japan).The Plecoptera (stonefly) data were sourced from Zhao et al. (2020).Inserted figure: Base substitution rate versus age diagram.This figure illustrates the base substitution rate over time.The thick red curve represents a power trendline, accompanied by its equation and an R 2 value of 0.0266.The red thin curve depicts the expected approximate curve around an estimated age of 320 Ma (Carboniferous).
TA B L E 2 Odonata species and the corresponding accession numbers.In FigTree v1.4.4,we can obtain outputs for the posterior age ("Node ages"), posterior probability ("posterior"), and "rate median" (mean of three rates for three branches at a specific node;

Species
not constant due to the application of a relaxed clock model; Drummond et al., 2006).These values are displayed at each node in Figure 1.Accordingly, we created a diagram (Figure 1 inset) illustrating the base substitution rate ("rate median," displayed at each node) versus age ("Node age," displayed at each node).This diagram includes a power trendline, its equation, and an R 2 value of 0.0266.
To examine a time-dependent base substitution rate influenced by the Quaternary calibration (Ho et al., 2005;Papadopoulo et al., 2010), we constructed a dated tree solely calibrated by the Quaternary dates of Q1-Q6 and J (Figure 4).We also generated a tree solely calibrated by the pre-Quaternary dates of A to I, K, and X, excluding Q1-Q6 and J (Figure 5).

| Dated tree
The posterior probabilities at the basal or crown nodes of Archaeognatha, Zygentoma, Palaeoptera, Ephemeroptera, Odonata, Epiprocta, Anisoptera, Zygoptera, and Plecoptera are all 1, except for one instance where the probability is 0.99.This high probability indicates the reliability of the current topology (Figure 1).
The crown age of Insecta was estimated to be 393.39Ma, with Archaeognatha (Apterygota) identified as the oldest lineage within TA B L E 3 Fossil and geological event calibrations.A multifurcation is observed among the species inhabiting the Rykyu-Taiwan-Japan Island region, including Chinese species.This complex pattern, reflecting a history of vicariance, began around 1.55 Ma and is observed in four ingroup species for Anisoptera (Q1-Q4) and two ingroup species for Zygoptera (Q5 and Q6; see Figure 3).The differentiation of Matrona and Calopteryx species is dated to 1.76 Ma, as calibrated by the point "J."

| Base substitution rates
The base substitution rate versus age diagram reveals that the rate is not constant, as seen in the application of the relaxed clock model.
Instead, it exhibits an exponential increase since approximately 20 Ma, as depicted in the equation shown in the inset (Figure 1   F I G U R E 4 BI Tree Calibrated Solely by Pre-Quaternary Fossil Dates, Utilizing BEAST v1.10.4.Observe that point K is indeed a fossil calibration point, but its date falls within the Quaternary.Both point K and points Q1-Q6 have the time to most recent common ancestor (tMRCA) set to default and have not been assigned specific dates.Notably, this tree does not incorporate any Quaternary node ages, with the minimum node age being 5.05 Ma, which closely aligns with the maximum root age of 5.8 Ma in Figure 5. Inset: Base substitution rate versus age diagram.It is essential to note that the rates depicted here are significantly slower than those shown in Figure 1, often falling into a single-digit range.
F I G U R E 5 BI tree calibrated solely by Quaternary dates, employing BEAST v1.10.4.It is important to recognize that point K is a fossil calibration point, and its date falls within the Quaternary, making it part of the current calibration.Calibration points A to I have the time to most recent common ancestor (tMRCA) set to default values, without specific dates assigned.Remarkably, this tree lacks the extended timescales seen in Figure 1, with the root age estimated at only 5.8 Ma, in stark contrast to the 393.39 Ma age in Figure 1.Inset: Base substitution rate versus age diagram.Take note that the rates displayed here are notably higher than those in Figures 1 and 4.
Figure 2 demonstrates that even mitochondrial genes with rapid base substitution rates, as discussed in Osozawa, Sato, et al. (2017), do not exhibit saturation toward ancient time periods.

| Impact of the Quaternary calibration on dated tree
The haplotype networks shown in Figure 3 support the ingroup or clade setting for the Quaternary calibration.
The phylogenetic tree calibrated exclusively using pre-Quaternary dates is depicted in Figure 4, and it notably exhibits a topology that aligns with the fully calibrated tree displayed in Figure 1.
However, the estimated ages for nodes in Figure 4

| Dated tree
It is important to note that our application of mostly mitochondrial genes successfully estimated ancient divergence times, and there is no evidence of mutation saturation (Figure 2).Our results support the paraphyly of Apterygota and the monophyly of Archaeognatha, Zygentoma, Pterygota, Palaeoptera, Ephemeroptera, Odonata, Epiprocta, Zygoptera, and each Odonata family.Ephemeroptera is sister to Odonata, Epiprocta is sister to Zygoptera, and Palaeoptera is sister to Neoptera (Figure 1; c.f., Figures 3 and 4; Table 1).These We presented a timetree for the East Asian Coeliccia damselfly, with details on vicariance provided in Osozawa, Sato, et al. (2017).
This timetree was calibrated exclusively using the Quaternary date of geological events, set at 1.55 ± 0.15 Ma.However, a significant issue emerged during this analysis, as the calculated basal node age (root age) was found to be much younger than what would be considered reasonable.This chronological problem is not unique to the Coeliccia damselfly but also appeared in the analyses of other East Asian Odonata timetrees, such as the Anotogaster dragonfly (Osozawa et al., 2013) and the Chlorogomphus dragonfly (Osozawa & Wakabayashi, 2015), as well as in other East Asian insects like the Papilio butterflies (Osozawa et al., 2013), Mycalesis butterfly (Osozawa, Takáhashi, et al., 2015), Pyrocoelia firefly (Osozawa, Oba, et al., 2015), Cicindela tiger beetle (Osozawa, Fukuda, et al., 2016;Osozawa, Ogino, et al., 2016), Ypthima butterfly (Osozawa, Takáhashi, et al., 2017), Platypleura cicada (Osozawa, Shiyake, et al., 2017), Ryukyu endemic five cicada group (Osozawa, Kanai, et al., 2021), and carabid ground beetle (Osozawa, Ogino, et al., 2016).Note that while these papers successfully describe extensive vicariance resulting from the isolation of islands from the Chinese continent, it has become evident that the inclusion of pre-Quaternary calibration points, in addition to Quaternary data, is necessary for achieving more precise dating in these analyses.
Epiophlebia superstes is a constituent of the Epiophlebioptera clade, with the sister group being the Anisopteromorpha, represented here by the Liassophlebiidae, a member of the stem group of the Anisoptera clade.Both of these groups collectively constitute the extant representatives of the Epiprocta major clade.Remarkably, E. superstes has persevered as an independent lineage, with an age as ancient as 205.75 Ma.Epiophlebia spp.have been identified in isolated regions encompassing the Japanese islands, China (the same species has also been reported in North Korea; Gunther et al., 2013), and the Himalayas.Despite this broad distribution, indications suggest that the vicariance in these areas has been relatively mild, implying that the isolation occurred recently and is not directly related to the age of this lineage (Busse et al., 2012).A parallel situation can be observed in the Araucariaceae conifer, an ancient clade with very young taxa in recent times (Escapa & Catalano, 2013).
Petaluridae, encompassing Tanypteryx pryeri and Petalura gigantea, was originally dated at 152.3 Ma, placing them in the earliest Cretaceous period.However, Ware et al. (2014) proposed that the entire Petalurida group, which includes Petaluridae, is an ancient clade with roots dating back to the Triassic (>201.3Ma).This extended timeframe for their diversification may be linked to the breakup of the supercontinent Pangaea.As detailed in our methods section (Appendix S2), it is noteworthy that the Pangaea breakup and the initiation of the Atlantic Ocean were recorded in the Santana Formation, occurring between 125 and 100 Ma (Martill, 2007; Table 3).This timeframe is notably more recent than the previously estimated crown age of 152.3 Ma for Petaluridae.Moreover, it is essential to acknowledge that the history of Pangaea's breakup is considerably more intricate, as we have explored in the context of mammalian evolution (Osozawa, 2023).
At the family level, differentiation within the extant Anisoptera occurred between 161.94 and 134.38 Ma, while for the extant Zygoptera, it took place from 116.4 to 3.85 Ma (Figure 1).Species-

| Quaternry vicariance increased biodiversity
Calibration points Q1 to Q6, representing Stylogomphus, Asiagomphus, Anotogaster, Chlorogomphus (Osozawa & Wakabayashi, 2015), Rhipidolestes, and Coeliccia (Osozawa, Sato, et al., 2017), were set at the 1.55 ± 0.15 Ma geologic event.Each of these multi-furcations or polytomies in the phylogenetic tree (Figure 1) and each of the haplotype network patterns (Figure 3) signify the emergence of multiple endemic species.The isolation of islands resulting from the opening of the Okinawa trough is a physical process, but it is evident that this process has significantly contributed to an increase in biodiversity.
While such isolation may lead to severe bottlenecks for endemic species, it is intriguing that this bottleneck effect does not align with the observed low genetic diversity and low nucleotide substitution rates found in this case, which contrasts with Zhai et al. (2017).
Calibration point J, representing Calopteryx and Matrona, was calibrated based on a fossil date of 1.76 ± 0.22 Ma.The divergence between Calopteryx japonica (found in Japan) and Matrona japonica (located in Amami-Okinawa; although Matrona basilaris in the Taiwan specimen could not be amplified, data from China were applied) appears to be a consequence or possibly a precursor to the 1.55 Ma vicariance event, as depicted in Figure 1.
In general, the process of back arc spreading, leading to the formation of continental islands, has the potential to trigger vicariance events that significantly contribute to increased biodiversity.There are numerous back arc basins in the western Pacific Ocean, some of which remain active, including the Okinawa trough.The notable high diversity of terrestrial organisms observed on continental islands, separated from continental landmasses due to sea-floor spreading of this nature, is likely a direct result of the physical isolation of these islands stemming from the rifting process, and the subsequent vicariance that ensues.
Conversely, oceanic islands are also created through back arc spreading, where volcanic edifices emerge in the ocean.However, these islands do not initially host terrestrial life because they surface above sea level as a result of ongoing volcanic activity.In such cases, the islands acquire their initial terrestrial species through dispersal from other landmasses.The species diversity associated with oceanic islands is more likely to increase through vicariance after the initial colonization via dispersal (Osozawa, Ito, et al., 2021;Osozawa, Kanai, et al., 2021;Osozawa, Ogino, et al., 2016).
We demonstrated that lotic damselflies tend to undergo vicariant speciation, in contrast to their lentic counterparts, as described in Osozawa, Sato, et al. (2017).A similar finding was reported by Letsch et al. (2016).However, it is important to note that within the Octogomphinae subfamily, there is no such tendency observed between the lotic Davidius and lentic Trigomphus.
Libellulidae, however, is predominantly composed of lentic species, where vicariance may be minimal, indicated by the similarity in genetic sequences between island populations, resulting in lower species diversity.An exception to this pattern is Sympetrum pedemontanum, which is a lotic species within the Sympetrum clade of Libellulidae.This lotic tendency is also observed in some species of mostly lentic Coenagrionidae damselflies.

| Increasing base substitution rate and biodiversity
Applying only the pre-Quaternary calibration resulted in the generation of very slow base substitution rates, as illustrated in the inset of  3).However, it is also important to consider an alternative perspective.The greater apparent biodiversity in geologically recent times might be influenced by the greater availability and preservation of, as well as the consequently increased number of studies conducted on, recent geological sections (as discussed in Rohde & Muller, 2005).Molecular phylogenetic analyses, in conjunction with paleobiological studies (Buatois & Mángano, 2018;Neige, 2015;Sahney et al., 2010), can be employed to test these alternative hypotheses (c.f., Barrier et al., 2001;Jønsson et al., 2012;Nakamura et al., 2021).
A possible driving factor behind the exponential increase in evolutionary rates and the likely rise in biodiversity toward the present time is extensive adaptive radiation (c.f., Ho et al., 2011).This could be attributed to the onset of glacial and interglacial cycles, coupled with environmental changes during the Quaternary period starting approximately 2.58 Ma, as discussed in Osozawa (2023).Following Ho et al. (2011), human mitochondrial DNA reflects adaptive changes in response to climatic variations (Mishmar et al., 2003;Ruiz-Pesini et al., 2004).Adaptive radiation within Murinae has been partly associated with positive selection and genomic changes (Roycroft et al., 2021).Genome analyses also suggest that the adaptive radiation of Heliconiini butterflies, the same target of Brower (1994), led to both phenotypic and genotypic variance (Cicconardi et al., 2023; see their Figure 1).Note that the Quaternary glaciations might have been triggered by the expansion of C4 land grasses and the evolution of sea diatoms, as this process led to increased carbon fixation, subsequently resulting in a decrease in atmospheric CO 2 concentration, as proposed by Taira (2007).The expansion of C4 grasses was a global phenomenon, encompassing regions such as North America and South America, and it initiated in the Oligocene and extended into the late Miocene, persisting to the present day.This expansion has relevance to the current glacial-interglacial period, although with some time lag, as highlighted by Cerling et al. (1997).
Although the available data plots were somewhat limited, we observed another peak in the base substitution rate during the Carboniferous to Permian period, as indicated in the inset of

ACK N OWLED G M ENTS
Chris Foote (Senior Editor), anonymous associate editor and tree reviewers contributed to enhancing the quality of the manuscripts.
We would like to extend our gratitude to Koji Tojo, Takuya Murata, Kunihiko Matsuhira, Yoshinori Kubota, Hidetoshi Sugita, Fumiyasu Sato, Takehiko Yamanaka, Junko Kobayashi, Kyoji Osozawa, and Akira Mishima for providing the analyzed samples that were invaluable for this study.Our special thanks go to Hyun-Yong Chung for Osozawa et al. (2012) illustrated how the back-arc spreading of the Okinawa Trough led to the formation of the Ryukyu Islands, resulting in their separation from the Chinese mainland.This seafloor spreading process initiated the islands' separation approximately 1.55 ± 0.15 Ma, and this isolation has persisted over time.
, while the Dominican amber was assigned to the Aquitanian period (ranging from 23.03 to 20.44 Ma; 21.735 ± 1.295 Ma; Iturralde-Vinent & MacPhee, 1996).They emphasized that "biostratigraphic information has been updated and adapted to the current geological time scale," following the most recent International Chronostratigraphic Chart established by the International Commission on Stratigraphy (Cohen et al., 2013; using v 2023/09 updates).Kohli et al. (2016) conducted an extensive examination of numerous crown dates available for Odonata fossil calibration (as presented in their Table figure.Therefore, these applications are not suitable for our current research objectives. E 2 Number of base changes in transition and tansversion versus corrected pairwide distance diagram for mitochondrial COI, COII, and 16S rRNA, as well as nuclear 28S rRNA.In MEGA 11, we considered the first and second codon positions and excluded the third codon position following the methodology ofYuan et al. (2022; MrBayesv3.2.6;Ronquist et al., 2012) and referred to the approach outlined inBi et al. (2023).

F I G U R E 3
Haplotype networks for Platycnemididae (Q6) and Rhipidolestidae (Q5).Haplotypes on the background half-moon represent a clade of Q6 and Q5, respectively, indicating vicariance associated with adaptive radiation.Dispersal is not considered.Refer to Figure 1 for corresponding details.Since our timetree covers a time span extending from 393.39 million years post-Silurian, calibrated up to 400.45 million years at point A, the base substitution rate versus age diagram captures extremely ancient variations in substitution rates.Another noteworthy observation is a mild peak in the rate between the Carboniferous and Permian periods, with rates of 0.0553 s/s/my at the crown node of Paleoptera (313.7 Ma) and 0.0412 s/s/my at the crown node of Odonata (284.24Ma).
Figure 5, the estimated node ages are notably younger, with the root age being only 5.8 Ma (resulting in every node age being less than 5.8 Ma).The base substitution rates in Figure 5 are substantially higher, reaching up to 7.56 s/s/my.Similar to Figure 4, there are two instances where branches in this tree appear to be collapsed or compressed, indicating inaccuracies in Figure 5.
Figure1, which incorporates both Quaternary and pre-Quaternary calibrations, is more reliable.Epiprocta is established as the sister to Zygoptera, while Epiophlebioptera is identified as the sister to Anisoptera.Within Anisoptera, Aeshnidae emerges as the sister to the remaining Anisoptera.Moving further, Petaluridae is recognized as the sister to Gomphidae, and the combined clade of Petaluridae + Gomphidae acts as the sister to the remaining Anisoptera.Additionally, Cordulegastridae is determined to be the sister to Chlorogomphoidea, and Cordulegastridae + Chlorogomphoidea is established as the sister to the remaining Anisoptera.Within this remaining group, Macromiidae + Corduliidae + Synthemistidae is the sister to Synthemistidae.These relationships are consistent withKohli et al. (2021) andSuvorov et al. (2022).These studies tend to overestimate node ages, and this could be attributed to factors such as the utilization of maximum age and the reliance on calibrations solely based on pre-Quaternary data.
level differentiation in Odonata occurred after 34.92 Ma, mainly during the Paleogene and Neogene periods.The endemic species and subspecies, which were calibrated with points Q1 to Q6 at 1.55 Ma, differentiated during the Quaternary.It is important to note that in the dated trees created by Kohli et al. (2021) and Suvorov et al. (2022), species-level differentiation is predominantly from the Paleogene, with some instances extending into the Cretaceous, and these estimates are now known to be overestimated due to the use of maximum age and reliance on calibrations solely based on pre-Quaternary data.Yu and Bu (2011) conducted a cladistic analysis of damselflies and demonstrated that Mesopodagrion and Rhipidolestes belong to different families, effectively splitting the former "Megapodagrionidae" into two distinct families, Megapodagrionidae and Rhipidolestidae.The two distinct clades of "Megapodagrionidae" in Figure 1 are indicative of this subdivision.

Figure 4 .
Figure 4. Conversely, when solely the Quaternary calibration was applied, it led to the emergence of very rapid base substitution rates,

Figure 1 .
Figure1.This peak may be analogous to the Late Paleozoic Ice Age, as discussed byMontañez et al. (2016), Rolland et al. (2019), andRosa and Isbell (2020).While glacial episodes in Earth's history are known to be influenced by factors such as landmass configuration, including the Gondwana supercontinent, it is possible that feedback from biological developments also played a role in initiating this glaciation.The Late Paleozoic glaciation is believed to have been triggered by the proliferation of terrestrial plants, specifically ferns, leading to the formation of thick coal layers during the Carboniferous period, as proposed byFranks et al. (2014).This process had the effect of increasing carbon fixation, which in turn effectively reduced atmospheric CO 2 levels, as previously suggested byTaira (2007) and corroborated byMontañez et al. (2016) andRolland et al. (2019).

| Geological evaluation of fossil dates for BEAST v1.10.4 analyses
thought to have acted as barriers to migration, thus triggering vicariance.